Molecular Modeling for Calculation of Mechanical Properties of 


Epoxies with Moisture Ingress 


T.C. Clancy*, S.J.V. Frankland* 
National Institute of Aerospace 
100 Exploration Way 
Hampton, 1 7 A 23666 


J.A. Hinkley^ T.S. Gates^ 
NASA Langley Research Center 
Hampton, VA 23681 


Abstract 

Atomistic models of epoxy structures were built in order to assess the effect of 
crosslink degree, moisture content and temperature on the calculated properties of a 
typical representative generic epoxy. Each atomistic model had approximately 7000 
atoms and was contained within a periodic boundary condition cell with edge lengths of 
about 4 nm. Four atomistic models were built with a range of crosslink degree and 
moisture content. Each of these structures was simulated at three temperatures: 300 K, 
350 K, and 400 K. Elastic constants were calculated for these structures by monitoring 
the stress tensor as a function of applied strain deformations to the periodic boundary 
conditions. The mechanical properties showed reasonably consistent behavior with 
respect to these parameters. The moduli decreased with decreasing crosslink degree with 
increasing temperature. The moduli generally decreased with increasing moisture 
content, although this effect was not as consistent as that seen for temperature and 
crosslink degree. 

1. Introduction 

Epoxies are useful compounds in many applications. In aerospace applications, 
they are used as adhesives as well as matrix components in composites. The effects of 
hygrothermal aging on epoxy materials are difficult to predict. Hygrothermal aging 
typically refers to the deterioration of the perfonnance of a material due to prolonged 
exposure to moisture and variable temperature conditions. In this paper, structural 
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features associated with hygrothermal aging are studied with atomistic level computer 
simulation. The assumption is made that variables such as crosslinking degree and 
moisture content will vary with time and exposure to enviromnental conditions. 

Water absorption in epoxies can be a complex process even under controlled 
laboratory conditions. Experimental studies have found water to exist in both “free” and 
“bound” states[l]. Additional studies have even found two different “bound” states for 
water in epoxy [2], Absorption of water can bring about both reversible and irreversible 
changes to the epoxy resin[3]. However, it is reasonable to assume that water content 
will generally increase with exposure to humid or wet enviromnents and reach a 
saturation value. Under certain conditions, the modulus of the resin can increase with 
absorbed water due to strong dipole-dipole interactions [2, 4]. Generally, however, water 
plasticizes the epoxy, reducing properties such as modulus and glass transition 
temperature, T g [5,6]. Chemical bonds can break due to hydrolysis and thermo-oxidative 
aging, effectively reducing the crosslinking degree in the cured epoxy [7]. Crosslinking 
reactions can also continue in an epoxy over time resulting in increased stiffness[8]. 

Given these complex relationships between the aging process and these two 
structural features of interest (moisture content and degree of crosslinking), the effect of 
hygrothermal aging is not directly studied. Rather, in this initial work, quantitative 
values for mechanical properties are predicted as a function of these structural features 
and temperature in the atomistic models. Although time dependence is not considered 
directly here, these molecular structural features known to vary with time and exposure 
are considered. Therefore, aging effects on these epoxy models are studied indirectly. 
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Fully atomistic simulations of macromolecular systems can be difficult to 
construct because of the range of length scales involved and the complexity of the 
bonding. Most fully atomistic simulations of polymers involve amorphous linear 
polymers and the chain lengths are typically oligomeric in scale. Crosslinked polymers, 
such as epoxies, have highly connected complex molecular architectures. Due to the 
complex structure of crossli nk ed polymers, only a limited amount of research has been 
conducted in the fully atomistic modeling of crossli nk ed macromolecular structures. A 
method was developed recently for building crosslinked atomistic structures[9]. This 
method involves the instantaneous bond formation of crossli nk s between reactant species. 
This basic method has been used with some modification for most atomistic simulations 
of crosslinked macromolecules[10-12]. 

The atomistic modeling of polymers for calculating mechanical properties has 
been researched for several decades. A static method was developed early on using 
molecular mechanics (MM) techniques[13]. Stress-strain simulations in molecular 
dynamics (MD) have been employed for some time[14]. However, perfonning these 
simulations is not a trivial task. It can be difficult to obtain mechanical properties from 
atomistic simulation that compare well with experimentally measured values due to the 
difficulty of preparing the models and the various parameters involved[15, 16]. A 
comparably high rate of strain is necessary in these simulations due to the short time 
scales of MD simulations as compared to experimentally accessible time scales. 

Atomistic MD simulations typically employ the integration of time steps 1 fs in duration. 
A feasible (computationally affordable) run might employ I 0 5 to perhaps 10 6 time steps, 
which would be a duration of 10"- 10' ps. 
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The organization of this paper is as follows. In section 2 of this paper, the 
methodology for building these atomistic models and calculating the properties is 
described. In section 3, results obtained from these simulations are presented. In section 4, 
the conclusions are presented. 

2. Methodology 

In this section, the methods used for generating the models and calculating the 
properties of these models are described. The details of the molecular species and the 
force field are listed. In section 2.1, the methodology for constructing the crosslinked 
network is described. In section 2.2, the method for calculating the elastic constants is 
presented. 

2.1. Modeling of epoxies and crosslinked polymers 

Fig. 1 shows the reactant molecules. The molecule containing the epoxy groups 
is the diglycidyl ether of bisphenol A (DGEBA) shown in Fig. la. There are two isomers 
for the reactant containing the crosslinking amine groups, diethyltoluene diamine 
(DETD A) shown in Figs. lb-c. These two isomers are nearly identical. The second 
amine group is bonded to the number 4 toluene ring carbon atom in the isomer in Fig. lb, 
while the second amine is bonded to the number 6 toluene carbon atom in the isomer in 
Fig. lc. There were 108 DGEBA molecules and 54 DETDA molecules in this simulated 
mixture for a total of 6966 atoms. There were 42 molecules of the isomer shown in Fig. 
lb and 12 molecules of the isomer shown in Fig. lc. This ratio comes from the known 
ratio of these isomers in commercially available samples of DETDA[17]. Fig. 2 shows 
the chemical reaction which takes place between the epoxide functional groups in the 
DGEBA reactant and the amine functional groups in the DETDA reactant. In Fig. 2a, a 


4 



primary amine group reacts with an epoxide group to form a chemical bond between the 
nitrogen of the amine and the terminal carbon of the epoxide group. The carbon-oxygen 
bond breaks between the terminal carbon and the epoxide oxygen, leaving an alcohol 
functional group. In Fig. 2b, the product of the reaction depicted in Fig. 2a undergoes 
further reaction with an additional epoxide group. The nitrogen atom now forms the 
center of a crosslinked molecular architectural structure. 

The AMBER force field[ 1 8, 19] was used in simulating these structures. The 
TIP3P[20] parameters were used to simulate the water in these simulations. Charges 
were estimated by using a bond decrement method used for the CVFF force field in the 
Materials Studio software[21]. This bond decrement method assigns equal magnitude, 
opposite sign charges to each pair of atoms connected by a bond. This produces a 
consistent, if approximate charge assignment. The LAMMPS[22] software was used to 
run these molecular mechanics (MM) and molecular dynamics (MD) simulations. 
Molecular mechanics typically refers to energy minimization techniques applied to static 
structures. Molecular dynamics involves the integration of the classical equations of 
motion of a system of atoms undergoing dynamical motions such as vibrations and 
rotations of atoms. A 1 fs time step was used in all the MD simulations. During 
equilibration runs, the LAMMPS software restricts the periodic boundary conditions to a 
rectilinear shape. The xyz dimensions can be coupled when a barostat is applied, which 
can be used to retain a cubic geometry. The simulation of the unreacted monomers was 
run with MD using a constant pressure-constant temperature (NPT) algorithm. The xyz 
box dimensions were coupled in order to maintain a cubic geometry. This was run for 
2x1 Cf time steps (200 ps). The final configuration was taken and used for the 
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construction of the crosslinked epoxy. Fig. 3 shows the atoms contained within the 
constraints of a periodic boundary condition cell. The edge length of this cubic cell is 
approximately 4 nm. The x and z coordinates lie in the plane of the figure. 

The final configuration of atomic coordinates of the reactant mixture is used to 
form a crossli nk ed network by the addition of chemical bonds consistent with the 
chemical reactions depicted in Fig. 2. Epoxide groups react with amine groups to form 
chemical bonds resulting in a complex network structure. In order to fonn the network, 
the distances between atoms which will be connected by crossli nk ed chemical bonds 
must be calculated. The first step is to calculate all the distances between all possible 
pairs of the terminal carbon atoms in the epoxide functional groups and the nitrogen 
atoms in the amine groups. The second step is to tabulate a list of potential bonds 
between the epoxy terminal carbons and the amine nitrogens, prioritizing by shortest 
distance. The bonding is restricted such that each epoxy terminal carbon is bonded to a 
maximum of one amine nitrogen and that each amine nitrogen is bonded to a maximum 
of two epoxy terminal carbons. The formation of these bonds will be restricted on the 
basis of a cutoff radius, r cuto ff, between these atom pairs. An r cuto ffX alue greater than 7-8 
A will typically result in configurations which are difficult to equilibrate with 
conventional MM and MD techniques. The final step is to form these bonds. The 
hydrogen atom that is attached to the reacting nitrogen is then bonded to the epoxide 
group, forming an alcohol as depicted in Fig. 2a. This new crosslinked structure is then 
used for the MM/MD simulation. 

The newly created atomistic configuration must be carefully equilibrated since the 
use of a relatively large cutoff radius, r cuto ff, results in unrealistic bond lengths. A short 
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initial energy minimization is applied to the structure using a conjugate gradient 
algorithm in the LAMMPS software with an energy delta tolerance of 10‘ 4 kcal/mol. 

This is followed by an NPT MD simulation at a temperature of 300 K and a pressure of 1 
atm for 2xl0 5 time steps (200 ps). The procedure described in the above paragraph can 
be repeated in order to increase the degree of crosslinking, a. The degree of crosslinking, 
a, is a unitless fraction with a maximum value of 1.0 for this stoichiometrically 
equivalent system. When a is equal to 1.0, all epoxide groups are bonded to an amine 
and each amine is bonded to 2 epoxide groups. By specifying r cuto jj and by choosing the 
number of iterations of the procedure, models with a range of a values can be obtained. 

Once the desired degree of crosslinking, a , is achieved, a more complete 
equilibration is applied. An energy minimization procedure is applied with a restriction 
of energy delta tolerance of 10" 4 kcal/mol. The structure was simulated with MD for 10 5 
time steps (10" ps) at 450 K using an NPT algorithm with no coupling between the three 
box dimensions. The dimensions remain nearly cubic however. Subsequently, the 
temperature was reduced linearly by 25 K over a simulation period of 20 ps. This was 
followed by 100 ps of MD simulation at 425 K. This procedure is repeated in stepwise 
manner down to 300 K. The final trajectory at each temperature is saved for running 
deformation simulations. 

In order to create atomistic models with moisture present, water molecules were 
added randomly to the structures. Water molecules were added at 5 and 9.5 wt. %. 
Following the addition of water molecules, energy minimization was applied. The 
structures were then subsequently equilibrated with MD from 450 K to 300 K using the 
stepwise cooling procedure described in the previous paragraph. 
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In order to obtain additional configurations for running deformation simulations, 
the MD simulation was run for an additional 100 ps. The final trajectory at each 
temperature is saved for running deformation simulations. The MD simulation was 
repeated for another 100 ps in order to obtain 3 separate atomistic configurations with 
their associated velocities for a given degree of crosslinking (a), moisture content (f wa t e , •), 
and temperature ( T ). Although 10 time steps ( 1 0“ ps) does not result in drastically 
distinct configurations, the properties calculated from these configurations did show 
statistical variation. 

One major statistical variation that was not included is the building of 
independent structures through the crosslinking reaction process. Since the reactant 
molecules were randomly mixed, the crosslinked structures could be built with different 
molecular configurations, despite having the same value of a. For the current work, 
however, with consideration of the computational resources required to run the MD 
simulations, multiple independent structures with identical parameter values were not 
built. 


2.2. Method for calculating elastic constants 

The method used here for calculating the elastic constants is similar to a method 
for static defonnation of periodic atomistic models[13], but applied to a dynamic 
simulation in a constant volume-constant temperature (NVT) ensemble. The 6x6 elastic 
constant matrix C is determined by the partial derivatives of the stress tensor, r, with 
respect to the defonnation, s, as indicated by Eq. 1 . 
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A series of deformations are applied to the periodic boundary cell in order to estimate 
values of these partial derivatives. These deformations are the following: A uniaxial 
deformation is applied to the x coordinate. The remaining two dimensions (yz) are 
unchanged. This is repeated for the y and z coordinates while keeping the remaining (xz 
and xy, respectively) coordinates fixed. Likewise, three shear deformations are applied. 
All deformations are in the positive direction, such that the periodic boundary coordinate 
x is increased during a uniaxial deformation along x for instance. These deformations are 
performed in an NVT ensemble using the LAMMPS simulation software. The 
deformation is applied in a continuous fashion at every time step at a constant rate with 
respect to the original dimensions. For most of the deformations, the periodic boundary 
cell was deformed by 5% over a simulation period of 10" ps (10 time steps). This is a 
relative deformation rate of 5xl0 8 s' 1 or 0.02 A/ps. Though high by experimental 
standards, these are typical rates of deformations for MD simulations [2 3, 24], The stress 
tensor was recorded at 10 time step (10 fs) intervals leading to a collection of 10 4 data 
points. This stress is calculated from the virial stress formula. These are usually 
averaged in batches of 10 time step intervals (10 ps). For example, in the case of the 
deformation of the x coordinate, the xx component of the stress tensor, r xx , is plotted as a 
function of the relative deformation, r c xx , and the slope of this line is obtained. This value 
is placed in the first column, first row in the C elastic constant matrix. From the same 
simulation, plots of z vy and r zr versus the degree of deformation, £ xx , are obtained. These 
are placed in the 2 nd and 3 ld row of the 1 st column of the C elastic matrix. The first three 
columns of the C matrix are therefore obtained from the three tension simulations in x, y 
and z. The last three columns are obtained from the three shear deformation simulations. 
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Therefore, the entire 6x6 C matrix is generated. For this amorphous polymer model, 
isotropic symmetry is assumed. A least squares fitting method is used to obtain the 
elastic constants, // and /t from Eqs. 2-6. [25] 


4a -2b + 3c 


n = 


33 

(2) 

2n + c-15// 

A = 


6 

(3) 

a — C n + C 22 + C 33 

(4) 

b — C 12 + C 13 + C 2I + C 23 + C 31 + C 32 

(5) 

c = C 44 + C 55 + C 66 

(6) 


The Young’s modulus, E, and shear modulus, G, can then calculated from Eqs. 7-8. 

3/1 + 2 u 
E = p — 

X + n (7^) 

G = <u (8) 

3. Results 

In this section, the results obtained from the methods applied in section 2 are 
presented. Details associated with the construction of the atomistic epoxy models are 
noted. This is followed by the results obtained from the applied defonnation of the 
equilibrated models. 

3.1. Construction of the epoxy models 

Four atomistic structures with different degrees of crosslinking, a, were built. 

The two highly crosslinked structures, (a= 0.86, 0.80) required two successive 
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crosslinking reactions as described in section 2.1. The two with relatively low degree of 
crosslinking (a= 0.56, 0.38) were obtained with only one crosslinking reaction. Fig. 4 
shows an example result of the degree of crosslinking, a, obtained by using various 
cutoff radius values, r cuto jf. The fdled squares indicate the value of a which would be 
obtained in the first crosslinking reaction as a function of r cuto ff. The open circles indicate 
the total value of a obtained as a function of r cuto ff after the first crosslinking was applied 
using a cutoff radius value of 8 A. This indicates that it is possible to obtain highly 
crosslinked structures with only two cycles of the procedure. 

The equilibration procedure described in section 2 was applied to all four epoxy 
models. Fig. 5 shows the density, p , calculated as a function of temperature, T, for one of 
these equilibration simulations ( a= 0.80, and f wate r = 0 wt. %). The density and 
temperature are sampled every 1 ps during the 100 ps run. The last 50 ps of these data 
values are averaged while the first 50 ps are ignored in order to allow the model to come 
to equilibrium at each new temperature. The error bars shown in Fig. 5 are obtained from 
the standard deviations of these averages. Although such plots are often used to estimate 
the glass transition temperature, T g , this is not attempted here. The results indicate the 
expected decrease in density with increasing temperature, with a final density at 300 K 
(1.08 g/cm ) slightly lower than typically seen experimentally (1.16 g/cm ) [26], The 
final density for the highest degree of crosslinking (a = 0.86, and f water = 0 wt %) is only 
slightly greater (1.09 g/cm ). Using the data from a restricted temperature range (T = 
300-375 K), the volumetric coefficients of thermal expansion (CTE) are calculated. The 
data are restricted in this way to avoid the glass transition temperature and the boiling 
point of water. Table 1 shows these values for the twelve structures that were simulated. 
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In the case of high degree of crosslinking and no moisture present (a = 0.86 ,f wate r = 0 wt. 
%), the CTE is 278 (xlO' 6 ) K" 1 . A recent experimental paper[17] lists a linear CTE of 
79.99 (xlO' 6 ) C' 1 for a DGEBA-DETDA epoxy. This would correspond to a volumetric 
CTE of 240 (xlO' 6 ) K' 1 which is in the same general range as the value in row one of 
Table 1. Generally, the CTEs increase with increasing moisture content and decreasing 
degree of crosslinking, although there are some inconsistencies. The most notable 
inconsistencies are with the (a = 0.80) values. The (a = 0.80 ,f wa ter = 0 wt. %) value is 
lower than that for (a= 0.86 ,f wate r = 0). Also, the CTE value actually decreases from (a 
= 0.80 , /water = 5 Wt. %) to (<2 = 0.80 ,f wa ter = 9.5 Wt. %) . 

Following equilibration of each of the four neat equilibrated epoxy structures, 
water molecules were added at both 5 wt. % (135 water molecules) and 9.5 wt. % (270 
water molecules). The equilibration procedure was repeated as described in section 2.2. 
A total of 12 distinct atomistic structures were created with the variation in a (= 0.86, 
0.80, 0.56, 0.38) and f wate r ( = 0, 5, 9.5 wt. %). Instantaneous snapshots of the positions 
and velocities of the atoms of the structures were preserved at a range of temperatures 
from the equilibration simulation. Snapshots were chosen from three temperatures (T = 
300 K, 350 K, 400 K) as starting points to run deformation simulations. 

3.2 Mechanical properties of the models 

In this section, results are presented for the calculation of elastic constants of the 
models. We expect to see the moduli decrease for increasing temperature ( T ), reduced 
degree of crosslink (a) and increased water content (f wa ter )• Generally, these trends are 
observed with some inconsistencies. As noted in section 1, the relationship between the 
aging process and these structural features ( a,f wate , ) can be quite complex even under 
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controlled laboratory conditions. Therefore, we do not discuss the effects of aging 
directly, but rather consider the results of the calculations as a function of the controlled 
structural variables. 

Fig. 6 shows a plot of the stress tensor component z xx as a function of the degree 
of strain deformation in the x direction, s xx . These values are averages over 1000 data 
points taken at 10 time step intervals. Therefore, each data point is averaged over a 
duration of 1 ps. These values over the 10 ps interval fall on a linear slope. The slope of 
this line is shown in Fig. 6, calculated from a least squares fit of the averaged values. 

The partial derivative (dr xx /ds xx ), is estimated as the slope of this line. 

A large number of relatively short MD simulations were run in order to collect the 
data for the mechanical properties. With 4 values of a, 3 values of f wa ter, and 3 values of 
T, there are 36 variations in the parameter set. With the 3 sample configurations for each 
parameter set, there are 108 distinct configurations for running deformation simulations. 
Each of these has six deformation simulations applied to it, for a total of 648 MD 
deformation simulations. As detailed in section 2.2, the results from the 6 deformations 
are averaged using a least squares procedure to obtain the two Lame constants, X and //, 
from which the elastic constants can be defined. The elastic constants from the three 
independent configurations with the same set of parameters (a, f wate r, T) are then 
averaged and the standard deviation is calculated. 

Fig. 7 shows the Young’s modulus, E, as a function of crosslink degree, a, 
calculated from atomistic models with no moisture content ( f water = 0 wt. %). The three 
lines are shown for data calculated at three temperatures: 300 K, 350 K and 400 K. The 
values shown are taken from the average of the three configuration runs for a given set of 
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parameter values. The error bars shown are taken from the standard deviations of these 
three data values. Some error bars are actually smaller than the symbols and are therefore 
not readily visable. Generally, there is an increase in E with increasing a, within the 
limits of the error bars. At the two highest values of a, there seems to be a decrease in E. 
This slight decrease is most likely due to statistical error. The maximal values for the 
two highest crosslinking densities at 300 K are 2.8 ± 0.2 GPa and 2.6 ±0.1 GPa. This 
compares well with recent experimental measurements of DGEBA-DETDA epoxy which 
was measured at 2.7GPa[27]. Fig. 8 shows the shear modulus, G, as a function of a for 
three different values of temperature. The trends follow those from Fig. 7, which is to be 
expected from the least squares fitting procedure described in section 2.2. 

Fig. 9 shows the Young’s modulus, E, as a function of the temperature, T, for 
three values of moisture content, f water = 0, 5 and 9.5 weight %. The crosslink density, a , 
is fixed at 0.80. The modulus declines with increasing T for all values of f wate r, however 
the trends are not completely as expected as the lines cross. Ideally, all values would 
decline consistently with increasing/,, ate r and decreasing a. It appears that it is much 
more difficult to get reliable or expected trends as a function of moisture content as 
compared to as a function of temperature and crosslink degree. This may be due to the 
limitations of the modeling capability and current hardware/software technology. 
However, experimental findings have suggested a complex structure of water within the 
epoxy[l, 2, 28]. This will require further analysis, but is beyond the scope of this paper. 

Fig. 10 shows the same plot as Fig. 9, but with a = 0.38. Here, the Young’s 
modulus decreases to very low values at 400 K, when moisture is present. For a = 0.80 
and no moisture, the modulus decreases by 35% with an increase in temperature from 
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300 K to 400 K, whereas the same change in conditions for a = 0.38 results in a 70% 
decrease in modulus. Considering the plot, it would appear that the low value of a is 
most influential in lowering E at higher temperatures as compared to the moisture content. 

Fig. 1 1 shows the Young’s modulus, E, as a function of moisture content ,f wa ter, 
for the epoxy simulations at 300 K. In most cases, a consistent decline in E with 
increasing water content is seen. The data for a = 0.86 are less consistent, however, 
which illustrates that these simulations are less reliable in calculating properties 
dependent on moisture content. Another important point is to note that the Young’s 
modulus, E, is greater for a = 0.80 than for a = 0.86 in all values of f wate r except 9.5 wt. 

%. This value at a = 0.86 and f water = 9.5 wt. % is the only case in this graph where E 
does not decrease consistently with increasing f wa ter- Fig. 12 shows the Young’s modulus, 
E, as a function of moisture content, f wate r, for the epoxy simulations at 350 K. In this 
case, the results appear to exhibit the expected dependencies (decline in E with 
decreasing a and increasing f wate ,) more closely than in Fig. 1 1 . Among the results that 
are not expected are the following: At f wate r = 0 wt. %, E is slightly lower for a = 0.80 
than for a = 0.86. For a = 0.80, E increases slightly as f mter increases from 5 to 9.5 wt. 

%, although, this increase is below the standard deviation calculated for either value. For 
a= 0.38, E increases slightly as f water increases from 0 to 5 wt. %, although, this increase 
is below the standard deviation calculated for either value. 

Fig. 13 shows the Young’s modulus, E, as a function of moisture content, f water , 
for the epoxy simulations at 350 K. A notable aberration is the slight increase in E as 
fwater increases from 0 to 5 wt. %, with a = 0.86. This increase is below the standard 
deviation calculated for E with f wa ter = 5 wt. % and a = 0.86. 
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4. Conclusions 


The motivation of this work was to study the difficult to predict effects of 
hygrothermal aging on epoxies. Because of the complicated nature of the relationship 
between these structural parameters ( oc,f wa t er ) and the aging process even under controlled 
laboratory conditions, the aging effects are studied indirectly. Atomistically detailed 
models of crosslinked epoxies were constructed with a range of degree of crosslinking ( a ) 
and moisture content (f water)- The elastic constants were calculated for these models. The 
moduli were found to decrease consistently with increasing temperature. There was a 
relatively consistent dependency of the moduli on the degree of crosslinking, a, with 
these moduli values decreasing with decreasing a. A less consistent dependency was 
found on the moduli as a function of moisture content. 

For each model with a specific value of parameters (. a,f wa ter , T), there were 3 
different samples. These samples were obtained as different configuration snapshots 
from a MD trajectory separated by I x 1 0 5 time steps (100 ps). Although this is a 
relatively short time span, the values obtained from these different configurations had 
statistical variation as indicated by the error bars in Figs. 7-13. Some additional 
averaging is obtained by the assumption of isotropic behavior and the fitting indicated by 
Eqs. 2-6. 

Although statistically independent configurations with identical parameter sets 
were not constructed, variation was addressed by the use of configurations generated 
from MD simulation trajectories. In addition, by assuming isotropic symmetry and 
averaging over deformations applied to the various coordinates, an additional structural 
averaging was achieved. 
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a 

f water (wt. /o) 

Volumetric 
CTE (K 1 ) 

0.86 

0 

278. 

5 

286. 

9.5 

360. 

0.80 

0 

188. 

5 

391. 

9.5 

375. 

0.56 

0 

360. 

5 

425. 

9.5 

477. 

0.38 

0 

405. 

5 

500. 

9.5 

584. 


Table 1. The volumetric coefficients of thermal, CTE, expansion. 
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Fig. 1. a) Diglycidyl ether of bisphenol A (DGEBA). b) 3, 5-Diethyltoluene-2, 4-diamine, 
c) 3, 5-Diethyltoluene -2, 6-diamine. (DETDA) 
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Fig. 2. The chemical reaction of an epoxide with an amine. 
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Fig. 3. The atoms of the epoxy simulation contained with the constraints of the periodic 
boundary conditions. The cartesian coordinate reference system is shown with respect to 
the structure. 
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Fig. 4. Degree of crosslinking, a, as a function of cutoff radius r cuto ff. 
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Fig. 5. The density, p , is plotted as a function of temperature, T, for a sample with a 

0.80, and f W ater = 0 wt. %. 



Fig. 6. The xx component of the stress tensor, r xx , as a function of the strain in the x 
direction, s xx . These data is taken from the a = 0.80 ,f wa ter = 0 wt. % system. 
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Fig. 7. The Young’s modulus, E, as a function of crosslink degree, a, for three different 
temperatures. These values are shown for samples with no moisture content (f water = 0 wt. 
%). 
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Fig. 8. The shear modulus, G, as a function of crosslink degree, a, for three different 
temperatures. These values are shown for samples with no moisture content (f wat er = 0 
%)• 



Fig. 9. The Young’s modulus, E, is plotted as a function of temperature, T, for three 
values of the moisture content (f wate r = 0, 5, 9.5 wt. %) and with degree of crosslinking, a 
= 0.80. 
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Fig. 10. The Young’s modulus, E, is plotted as a function of temperature, T, for three 
values of the moisture content (f wate r = 0, 5, 9.5 wt. %) and with degree of crosslinking, 
= 0.38. 



Fig. 11. The Young’s modulus, E, as a function of moisture content, f wate r, at 300 K. 
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Fig. 12. The Young’s modulus, E, as a function of moisture content, at 350 K. 
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Fig. 13. The Young’s modulus, E, as a function of moisture content, f wate r, at 400 K. 
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